Image.open("/home/xpwang/Downloads/colposcope3/train/benign/20130902133200.jpg")
mixture of Gaussian the EM is the following.
Repeat until converged:
E-step. Compute the expectations, or responsibilities $r_{j,i}$, of the latent variables $v_i$, where $j$ represents the classes, $y$ the data points, and $i$ the $ith$ points: $$\begin{equation} r_{j,i} = \frac{\omega_j P(y_i~|~\mu_j,\Sigma_j)}{\sum_j \omega_j P(y_i~|~\mu_j,\Sigma_j)} \end{equation} $$
M-step. Compute the maximum likelihood parameters given these responsibilities:
(1) The weigths for the mixing proportions, $\omega_j$, where $N$ is the number of data points: The estimate of the probability of a class as an outcome is the proportion of times it appears: $$\begin{equation} \omega_j ← \arg\max \sum_i \sum_j r_{j,i} log P(z_i = j) = \sum_i \sum_j r_{j,i} log \omega_j = \frac{1}{N}\sum_i r_{j,i} \end{equation} $$
(2) The means and the variances for the Gauss distributions, $\mu_j$, $\sigma_j$ : We want to maximize this term: $$\begin{equation} \arg\max \sum_i r_{j,i} log P(y_i~|~v_i = j) = \sum_i r_{j,i} log Gaussian(y_i; \mu_j, \sigma_j) \end{equation} $$ $$\begin{equation} \mu_j ← \frac{1}{\sum_i r_{j,i}} \sum_i r_{j,i} y_i \end{equation} $$ $$\begin{equation} \sigma_j ← \frac{1}{\sum_i r_{j,i}} \sum_i r_{j,i} (y_i - \mu_j)^2 \end{equation} $$

def multile_stage(image,I_lamda = 0.7,S_lamda = 0.1) :
# 转换成array
arr = np.array(image)
# 计算intensity 和 sataration
arr_I = np.mean(arr,axis =2)
arr_S = 1 - np.min(arr,axis=2) / (arr_I+1e-7)
# 根据阈值生成mask
ise = arr_I > I_lamda * np.max(arr_I)
sse = arr_S < S_lamda * np.max(arr_S)
# 计算梯度
magnitude = generic_gradient_magnitude(arr_I,sobel)
mse = magnitude > np.percentile(magnitude,97)
#第一阶段的效果
select_1 = ise & sse
print(ise.sum(),sse.sum(),mse.sum())
arr_1 = arr.copy()
arr_1[:,:,0][select_1] = 0
arr_1[:,:,1][select_1] = 0
arr_1[:,:,2][select_1] = 0
# 第二阶段的效果
select = ise & sse | mse
arr[:,:,0][select] = 0
arr[:,:,1][select] = 0
arr[:,:,2][select] = 0
return arr_1,arr,select_1,select,mse
def gassin_process(image,mask):
arr = np.array(image)
mask = mask.copy()
arr_I = np.mean(arr,axis =2)
arr_S = 1 - np.min(arr,axis=2) / (arr_I+1e-7)
hsv = color.rgb2hsv(arr)
S = arr_S[mask]
V = hsv[:, :, 2][mask]
X = np.vstack((S, V)).T
model = GaussianMixture(n_components=2, covariance_type="full", max_iter=20000)
model.fit(X)
resoult = model.predict(X)
mask[mask][resoult == 0] = False
arr[:, :, 0][mask] = 0.0
arr[:, :, 1][mask] = 0.0
arr[:, :, 2][mask] = 0.0
return arr,mask
def fill_process(image,mask,kenel_size = 20) :
image = np.array(image)
mask = mask.copy()
mask = dilation(mask,selem=square(5))
image[:, :, 0][mask != 0] = 0.0
image[:, :, 1][mask != 0] = 0.0
image[:, :, 2][mask != 0] = 0.0
while(mask.sum() != 0) :
init = mask.sum()
shape = mask.shape
print(mask.sum())
for i in range(mask.shape[0]) :
for j in range(mask.shape[1]) :
nz = 0
sum_nz = np.array([0,0,0])
if mask[i,j] != 0 :
for ni in np.array(range(kenel_size)) - kenel_size/2:
x = int(i + ni)
for nj in np.array(range(kenel_size))- kenel_size/2 :
y = int(j + nj)
if x >0 and x < shape[0] and y > 0 and y < shape[1] and image[x,y,:].sum() != 0 and image[x,y,:].sum() < 245*3:
nz += 1
sum_nz += image[x,y,:]
if nz != 0 :
image[i,j,:] = sum_nz / nz
mask[i,j] = False
# print(mask.sum())
if init == mask.sum() :
# image[:, :, 0][self.select] = np.median(self.arr[:,:,0])
# self.arr[:, :, 1][self.select] = np.median(self.arr[:,:,1])
# self.arr[:, :, 2][self.select] = np.median(self.arr[:,:,2])
break
return image
plt.figure(figsize=(40,40))
plt.subplot(1,2,1)
plt.imshow(image)
plt.title("Image")
plt.subplot(1,2,2)
plt.imshow(test1)
plt.title("after Intensity and Saturation selection")
<matplotlib.text.Text at 0x7fe96b5d9c50>
plt.figure(figsize=(40,40))
plt.subplot(1,2,1)
plt.imshow(test2)
plt.title("after gradient magtitude selection")
plt.subplot(1,2,2)
plt.imshow(test3)
plt.title("after GMM")
<matplotlib.text.Text at 0x7fe96b3b49b0>
plt.figure(figsize=(10,10))
# plt.subplot(1,2,1)
plt.imshow(test4)
plt.title("after fill process")
<matplotlib.text.Text at 0x7fe969383c18>


contour_list = find_contours(modual_2(),0.1)
plt.figure()
fig,ax = plt.subplots(figsize=(10,10))
ax.imshow(modual_2(),cmap = "gray")
for contour in contour_list :
ax.plot(contour[:,1],contour[:,0])
<matplotlib.figure.Figure at 0x7fe94bdb3630>
image
def sr_mark(image) :
arr = np.array(image)
sele_var = np.var(arr,axis=2)
# sele_min = np.min(arr,axis = 2)
choose = sele_var.copy()
choose[:,:] = False
choose[50:120,80:570][sele_var[50:120,80:570] <100] = True
choose[120:,80:643][sele_var[120:,80:643] <100] = True
# choose[60:633,63:575][sele_min[60:633,63:575] > 180] = True
arr[:,:,0][choose != 0] = 0
arr[:,:,1][choose != 0] = 0
arr[:,:,2][choose != 0] = 0
# arr[:,:,0] = median_filter(arr[:,:,0],size=5)
# arr[:,:,1] = median_filter(arr[:,:,1],size=5)
# arr[:,:,2] = median_filter(arr[:,:,2],size=5)
choose = dilation(choose,selem=square(5))
return arr,choose
def fill_process(image,mask,kenel_size = 20) :
mask = mask.copy()
image = np.array(image)
image[:, :, 0][mask != 0] = 0.0
image[:, :, 1][mask !=0] = 0.0
image[:, :, 2][mask != 0] = 0.0
while(mask.sum() != 0) :
init = mask.sum()
shape = mask.shape
for i in range(mask.shape[0]) :
for j in range(mask.shape[1]) :
nz = 0
sum_nz = np.array([0,0,0])
if mask[i,j] != 0 :
for ni in np.array(range(kenel_size)) - kenel_size/2:
x = int(i + ni)
for nj in np.array(range(kenel_size))- kenel_size/2 :
y = int(j + nj)
if x >0 and x < shape[0] and y > 0 and y < shape[1] and image[x,y,:].sum() != 0 and image[x,y,:].sum() < 245*3:
nz += 1
sum_nz += image[x,y,:]
if nz != 0 :
image[i,j,:] = sum_nz / nz
mask[i,j] = False
# print(mask.sum())
if init == mask.sum() :
# image[:, :, 0][self.select] = np.median(self.arr[:,:,0])
# self.arr[:, :, 1][self.select] = np.median(self.arr[:,:,1])
# self.arr[:, :, 2][self.select] = np.median(self.arr[:,:,2])
break
return image
def background(image,back,percent = 99.5) :
back = back.copy()
image = image.copy()
back[:,:,0] = median_filter(back[:,:,0],size=10)
back[:,:,1] = median_filter(back[:,:,1],size =10)
back[:,:,2] = median_filter(back[:,:,2],size= 10)
mask = np.min(image / back ,axis =2 )
image[mask > np.percentile(mask,percent),:] = back[mask > np.percentile(mask,percent),:]
return image,mask,back
_,mask = sr_mark(image)
paint_2 = fill_process(image,mask)
paint_4,mask_4,back_4 = background(paint_2,paint_2)
plt.figure(figsize=(8,8))
plt.imshow(paint_2)
plt.title("after var remove and fill")
<matplotlib.image.AxesImage at 0x7effaf738ac8>
plt.figure(figsize=(8,8))
plt.imshow(paint_4)
plt.title("after background fillprocess")
<matplotlib.image.AxesImage at 0x7f0e1eb84320>
show(1)
show(2)
show(3)
show(4)
show(5)
show(9)
show(np.random.randint(1)%80)
show(np.random.randint(80))
show(np.random.randint(80))